// G4提供的头文件
#include "G4RunManager.hh"
#include "G4SystemOfUnits.hh"
#include "G4UnitsTable.hh"
#include "G4AccumulableManager.hh"
#include "G4Run.hh"
#include "G4LogicalVolumeStore.hh"
#include "G4LogicalVolume.hh"
#include "G4AnalysisManager.hh"
// 用户编写的头文件
#include "G4ImproveRunAction.hh"
#include "G4ImproveDetectorConstruction.hh"
#include "G4ImprovePrimaryGeneratorAction.hh"

// 构造函数
G4ImproveRunAction::G4ImproveRunAction(){
    
    // 定义新的剂量单位
    const G4double milligray = 1.e-3*gray;
    const G4double microgray = 1.e-6*gray;
    const G4double nanagray  = 1.e-9*gray;
    const G4double picogray  = 1.e-12*gray;

    // 将定义的新单位添加到G4单位制系统中
    new G4UnitDefinition("milligray","milliGy","Dose",milligray);
    new G4UnitDefinition("microgray","microGy","Dose",microgray);
    new G4UnitDefinition("nanagray" ,"nanoGy" ,"Dose",nanagray);
    new G4UnitDefinition("picogray" ,"picoGy" ,"Dose",picogray);

    // 创建能量沉积管理系统
    auto accumulableManager = G4AccumulableManager::Instance();
    accumulableManager->RegisterAccumulable(fEdep);     // 将储存沉积能量的变量添加到管理器中
    accumulableManager->RegisterAccumulable(fEdep2);    // 将储存沉积能量平方的变量添加到管理器中

    // 获取分析器管理实例
    auto analysisManager = G4AnalysisManager::Instance(); 
    // 设置输出文件格式为root      
    analysisManager->SetVerboseLevel(1);
    analysisManager->SetDefaultFileType("root");
    analysisManager->SetFileName("../out/Spectrum");

}
// 析构函数
G4ImproveRunAction::~G4ImproveRunAction(){}

// 每次run开始之前运行
void G4ImproveRunAction::BeginOfRunAction(const G4Run*){

    // 不保持随机数种子，每次运行都是一次新的模拟
    G4RunManager::GetRunManager()->SetRandomNumberStore(false);

    // 创建能量沉积的管理类
    auto accumulableManager = G4AccumulableManager::Instance();
    // 每次开始前重置fEdep和fEdep2的值
    accumulableManager->Reset();

    auto analysisManager = G4AnalysisManager::Instance();
    analysisManager->OpenFile();

    analysisManager->CreateNtuple("Tree","Detector Data");      // 创建一个新Ntuple
    analysisManager->CreateNtupleDColumn("Edep");       // 第一个分支
    analysisManager->CreateNtupleDColumn("x");          // 第二个分支
    analysisManager->CreateNtupleDColumn("y");          // 第三个分支
    analysisManager->CreateNtupleDColumn("z");          // 第四个分支
    analysisManager->FinishNtuple();            // 结束创建

    analysisManager->CreateH1("Edep","Energy deposition",1024,0,3*MeV);     // 创建一个直方图
    analysisManager->CreateH1("E","energy",256,0.1,1*MeV);                  // 创建第二个直方图


}

// 每次run结束后运行
void G4ImproveRunAction::EndOfRunAction(const G4Run* run){

    // 获取当前运行中的事件数量,如果没有事件，则直接返回
    G4int nofEvents = run->GetNumberOfEvent();
    if(nofEvents == 0)  return;

    // 创建一个能量沉积管理器，该管理器可以处理能量沉积的累积、重置等
    auto accumulable = G4AccumulableManager::Instance();
    // 合并多线程模拟的结果
    accumulable->Merge();

    G4double edep  = fEdep.GetValue();      // 定义一个记录运行期间总能量沉积的变量
    G4double edep2 = fEdep2.GetValue();     // 定义一个记录运行期间总能量沉积平方和的变量

    // 计算能量沉积的均方根值rms
    G4double rms = edep2 - edep * edep / nofEvents;
    if(rms > 0.0)
        rms = std::sqrt(rms);
    else
        rms = 0.0;

    // 从runManager中获取用户的建模信息
    const G4ImproveDetectorConstruction* detConstruction = 
        static_cast<const G4ImproveDetectorConstruction*>(G4RunManager::GetRunManager()->GetUserDetectorConstruction());
    // 从建模信息中，获取敏感体积与其质量   
    G4double mass = detConstruction->GetScoringVolume()->GetMass();

    // 通过得到的能量沉积edep与敏感体积的质量mass，计算得到总剂量dose
    G4double dose = edep / mass;
    // 通过前面计算得到的均方根植rms与质量mass，计算得到均方根剂量
    G4double rmsDose = rms / mass;

    G4String runCondition;      // 定义一个字符串，用于储存获得的粒子类型

    // 从runManager中获取初级事件配置信息
    const G4ImprovePrimaryGeneratorAction* generatorAction = 
        static_cast<const G4ImprovePrimaryGeneratorAction*>(G4RunManager::GetRunManager()->GetUserPrimaryGeneratorAction());

    if(generatorAction){

        // 从初级事件配置信息中，获取粒子发射系统
        const G4ParticleGun* gps = generatorAction->GetGPS();

        // 获取发射的粒子类型的名称并添加到描述字符串中
        runCondition += gps->GetParticleDefinition()->GetParticleName();
        runCondition += " of ";

        // 获取发射的粒子的能量，并以最佳单位形式添加到描述字符串中
        G4double particleEnergy = gps->GetParticleEnergy();
        runCondition += G4BestUnit(particleEnergy,"Energy");

    }

    auto analysisManager = G4AnalysisManager::Instance();
    // 写入并关闭root文件
    analysisManager->Write();
    analysisManager->CloseFile();

    /***********************打印输出结束标志**********************************/
    if (IsMaster()) {
        // 如果当前线程是主线程（Master），则输出全局运行结束的消息
    G4cout
     << G4endl
     << "--------------------End of Global Run-----------------------";
    }
    else {
        // 如果当前线程不是主线程（Master），则输出本地运行结束的消息
    G4cout
     << G4endl
     << "--------------------End of Local Run------------------------";
    }

    /***********************数据的打印输出**********************************/
    G4cout
        << G4endl
        << " The run consists of " << nofEvents << " "<< runCondition           // 输出本次运行了几个事件
        << G4endl
        << " Cumulated dose per run, in scoring volume : "
        << G4BestUnit(dose,"Dose") << " rms = " << G4BestUnit(rmsDose,"Dose")   // 输出沉积的能量值与均方根值
        << G4endl
        << "------------------------------------------------------------"
        << G4endl
        << G4endl;

}

// 能量累加函数
void G4ImproveRunAction::AddEdep(G4double edep){
    
    // 累加每个event沉积的能量
    fEdep += edep;
    // 累加每个event沉积能量的平方和
    fEdep2 += edep*edep;
}
